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ABSTRACT 

The white dwarf luminosity function is an important tool for the study of the so- 
lar neighborhood, since it allows the determination of the age of the Galactic disk. 
Over the years, several methods have been proposed to compute galaxy luminosity 
functions, from the most simple ones — counting sample objects inside a given vol- 
ume — to very sophisticated ones — like the C~ method, the STY method or the 
Choloniewski method, among others. However, only the 1/Vmax method is usually em- 
ployed in computing the white dwarf luminosity function and other methods have not 
been applied so far to the observational sample of spectroscopically identified white 
dwarfs — in sharp contrast with the situation when galaxy luminosity functions are 
derived from a large variety of samples. Moreover, the statistical significance of the 
white dwarf luminosity function has also received little attention and a thorough study 
still remains to be done. In this paper we study, using a controlled synthetic sample 
of white dwarfs generated using a Monte Carlo simulator, which is the statistical sig- 
nificance of the white dwarf luminosity function and which are the expected biases. 
We also present a comparison between different estimators for computing the white 
dwarf luminosity function. We find that for sample sizes large enough the 1/Vmax 
method provides a reliable characterization of the white dwarf luminosity function, 
provided that the input sample is selected carefully. Particularly, the 1 / Vmax method 
recovers well the position of the cut-off of the white dwarf luminosity function. How- 
ever, this method turns out to be less robust than the Choloniewski method when 
the possible incompletenesses of the sample are taken into account. We also find that 
the Choloniewski method performs better than the 1 /Vmax method in estimating the 
overall density of white dwarfs, but misses the exact location of the cut-off of the 
white dwarf luminosity function. 

Key words: stars: white dwarfs — stars: luminosity function, mass function — 
Galaxy: stellar content — methods: statistical 



1 INTRODUCTION 

The white dwarf luminosity function is perhaps one of the 
most useful tools for deriving important properties of the 
solar neighborhood. In particular, but not only, the disk 
white dwarf luminosity function carries valuable information 
about the age of the Galaxy (Winget et al. 1987; Garci'a- 
Berro et al. 1988; Hernanz et al. 1994; Richer et al. 2000) and 
of the stellar formation rate (Noh & Scalo 1990; Di'az-Pinto 
et al. 1994; Isern et al. 1995; Isern et al. 2001). Additionally, 
the luminosity function of white dwarfs provides an inde- 
pendent test of the theory of dense plasmas (Segretain et 
al. 1994; Isern et al. 1997). Finally, the white dwarf lumi- 



nosity function directly measures the current death rate of 
low- and intermediate-mass stars in the local disk. Conse- 
quently, a reliable determination of the observational white 
dwarf luminosity function is of the maximum interest. 

Previous observational efforts, like the Palomar Green 
Survey (Green et al. 1986) have provided us with an in- 
valuable wealth of good quality data. Moreover, ongoing 
projects like the Sloan Digital Sky Survey (York et al. 2000; 
Stoughton et al. 2002; Abazajian et al. 2003, 2004), the 2 Mi- 
cron AU Sky Survey (Skrutskie et al. 1997; Cutri et al. 2003), 
the SuperCosmos Sky Survey (Hambly et al. 2001a; Hambly, 
et al. 2001b; Hambly, Irwin & MacGillivray 2001), the 2dF 
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QSO Redshift Survey (Vennes et al. 2002), the SPY project 
(Pauli et al. 2003), and others will undoubtely increase the 
sample of spectroscopically identified white dwarfs with reli- 
able determinations of parallaxes and proper motions, which 
are essential for an accurate determination of the white 
dwarf luminosity function. Last but not least, future space 
missions like Gaia (Perryman et al. 2001) will increase the 
sample of known white dwarfs with very accurate astromet- 
ric determinations even further (Torres et al. 2005), thus 
allowing a precise and reliable determination of the proper- 
ties of the disk white dwarf population. 

Over the years, several methods have been used to de- 
termine luminosity functions for all sort of objects, ranging 
from main sequence stars to galaxies. These include the most 
simple ones (counting stars inside a given volume) to very so- 
phisticated ones — like the C~ method (Lynden-Bell 1971), 
the STY method (Sandage et al. 1979) and the Choloniewski 
method (Choloniewski 1986). In spite of the variety of meth- 
ods currently used to estimate galaxy luminosity functions, 
the 1/Vmax method (Schmidt 1968) is the most commonly 
used method for estimating white dwarf himinosity func- 
tions, though, to the best of our knowledge, nobody has yet 
assesed in depth its statistical reliability for such a purpose. 
More specifically, up to now only two works have studied 
how good the 1/Vmax method performs in estimating the 
white dwarf luminosity function. In particular, Wood & Os- 
walt (1998) demonstrated by using a Monte Carlo simu- 
lator that the 1/Vmax method for proper-motion selected 
samples is a good density estimator, although it shows im- 
portant statistical fluctuations when estimating the slope of 
the bright end of the white dwarf luminosity function. Later 
on, Garci'a-Berro et al. (1999), using another Monte Carlo 
simulator, corroborated the previous study and, moreover, 
showed that the standard procedure used by the 1/Vmax 
method to asign error bars severely understimates the size 
of the real error bars for a typical sample of 200 objects. 
Additionally these last authors also showed that there was 
a bias in the derived ages of the solar neighborhood, conse- 
quence of the binning procedure. However, the most appar- 
ent conclusion of both papers is that, in general, selection 
effects or. simply, the inherent characteristics of the sample 
under consideration have a strong effect on the shape of the 
estimated white dwarf luminosity function, despite using an 
unbiased estimator, like the 1/Vmax method. 

In this paper we assess the statistical significance of the 
observational white dwarf luminosity function. For such a 
purpose we will use a controlled synthetic sample of white 
dwarfs generated with our Monte Carlo simulator. We dis- 
cuss in depth which are the typical biases introduced by the 
procedures used to select the sample. This includes both the 
bias in retrieving the correct slope for the monotonically in- 
creasing branch of the white dwaxf luminosity function and, 
most importantly, the bias obtained when retrieving the pre- 
cise location of the observed cut-off. This last point is of 
the maximum interest, since the observed drop-off of the 
white dwarf luminosity is currently one of the best estimar 
tors used to date the local neighborhood. We also present 
an independent estimate of the size of the error bars. Fi- 
nally we discuss the advantages and shortcomings of the 
several methods that can be used to obtain the observa- 
tional white dwarf luminosity function and we present a set 
of recommendations. The reader should take into account 



that in the present paper we only discuss the sampling bi- 
ases and do not take into account the measurement errors. 
Clearly, the effects of the measurement errors will affect the 
sampling biases and viceversa. Moreover, the effects of the 
measurement errors could bo as important as the sampling 
biases, although this remains still to be studied. Such an 
study is under preparation and will be published elsewhere. 
The paper is organized as follows. In section 2 we describe 
the different estimators most commonly used today for ob- 
taining luminosity functions. Section 3 is devoted to describe 
the Monte Carlo simulations used to compare the different 
methods previously described in §2, whereas in Section 4 we 
apply the different estimators to our Monte Carlo samples. 
Finally in Section 5 we summarize our major findings and 
we draw our conclusions. 



2 THE MOST COMMONLY USED 

LUMINOSITY FUNCTION ESTIMATORS 

2.1 Schmidt's estimator for proper motion and 
magnitude selected samples 

This method, also known as the 1/Vmax method, was first 
introduced by Schmidt (1968) in the studies of the quasar 
population. Later on Schmidt (1975) extended it to proper 
motion selected samples and Foltcn (1976) made a gener- 
alization of the method introducing the dependence on the 
direction of the sample. This turns out to be useful when 
studying stellar samples because the scale height of the 
Galactic disk introduces some biases. 

Consider a sample of stars having a lower proper mo- 
tion limit ytii and faint apparent magnitude limit mi, the 
maximum distance for which an object can contribute to 
the sample is: 

W = min [7r-^(Ai/Mi);7r-'lO°-'('"f-'">] (1) 

where n is the stellar parallax, /i is the proper motion, and 
m the apparent magnitude. If the sample is only complete 
to a certain upper proper-motion limit fj,„ and to a bright 
magnitude limit ruh, then there is also a minimum distance 
for which an object contributes to the luminosity function: 

rmin = max [7r-^(M/Mu); Tf-ho"-'''"''-"')] (2) 

Additionally, if the sample only covers a fraction of the sky, 
/?, then the maximum volume in which a star can contribute 
is: 

Vmax — "^/?(^max ^ ^'min) (3) 

The luminosity function, c^m, is then built by binning 
the sample in i € (l,iV) magnitude bins and performing a 

weighted sum over the objects in each magnitude bin, Ni. 
The weight with which every object contributes to the sum 
depends on the maximum volume in which it could be de- 
tected given its apparent magnitude and proper motion: 

j — l *^rnax 

Though this estimator is based on heuristic apprecia- 
tions about how a good estimator should be, it has been 
shown that it is unbiased (Felten 1976). However, the fact 
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that the estimator is unbiased does not guarantee a good 
estimate of the luminosity function if the input data is not 
properly selected. More specifically, the input sample should 
be complete and, additionally, Takeuchi et al. (2000) have 
pointed out that the 1 /Vmax method for magnitude selected 
samples is seriously affected when the input data — even if 
complete — is clustered or, more generally speaking, when 
the input data are not homogenously scattered. This, in 
turn, affects the derived shape of the luminosity function. 
Consequently, the 1 /Vmax method should only be used when 
homogeneity and completitude of the sample under consid- 
eration are guaranteed. This, of course, is not an easy task, 
and for most of the observational samples is an "a priori" 
assumption. 



2.2 Maximum Likelihood estimators based on the 
probability of selection 

Maximum likelihood estimators are based on the probabil- 
ity of selecting a given object in a sample and have been 
shown to be insensitive to sample inliomogenities. Moreover, 
by definition, these estimators are unbiased and have min- 
imum variance for large samples. Two variants have been 
already developed. The first one is a parametric estimator 
(Sandage et al. 1979), hereafter the STY estimator, whereas 
the corresponding non-parametric version (Efstathiou et 
al. 1988) is called the StepWise Maximum Likelihood es- 
timator (SWML). Both estimators are designed for magni- 
tude selected samples and their performances when eval- 
uating galaxy luminosity functions have been thoroughly 
tested using detailed Monte Carlo simulations (Willmer 
1997; Takeuchi et al. 2000). 

Following Luri et al. (1996) we define the likelihood 
function, as the product of the probability distributions 
of the variables of interest, under the assumption that all of 
them are statistically independent: 



jC{e) = llM{Mk\e) 



(8) 



(5) 



k=l 



where a; is a random variable with probability density 
'D{xk\0) depending on a set of unknown parameters 6 = 

{61,62, ■■■ ,0n) and realizations {xi, . . . , Xn^)- The value of 
6 that maximizes this function is the maximum likelihood 
estimator, 0ml, of the parameters. 

Observational selection may be modelled by introduc- 
ing a new probability density, M{xk\0), with the help of a 
normalization constant, Cj(^/, which depends upon the max- 
imization variables, and of a selection function, S{xk)- 

M{xk\6) ^CM'V{xk\6)Sixk) (6) 

A typical example of such a selection function is that 
obtained for the case of a sample which is complete up to a 
certain limiting apparent magnitude, miim- In this case the 
selection function is simply a Heaviside function, S{mk) = 
G(mfc — miim). Writing down the probability distribution 
as the product of the densities of the variables of interest 
in our sample — absolute magnitude, M, parallax, vr, and 
tangential velocity, vtan — we get the structure of the STY 
and SWML estimators for magnitude selected samples: 



It becomes obvious from the definition of likelihood that 

the statistical independence of the variables makes the max- 
imization process insensitive to the distribution of veloci- 
ties and to the parallax probability density. Note, however, 
that for the case in which we have a mixture of populations 
(thin and thick disk and halo, for instance) the velocities 
and the absolute magnitude are no longer independent vari- 
ables. The SWML method is then obtained by adopting for 
as a stepwise function: 

JV 

(^(M) = 'f>iW{Mi - M) (9) 

i=l 

where the window function W{Mi — M), is defined by: 

W{M, - M) = i I . 2 2 (10) 

1^ otherwise 

and (fii yields the luminosity function of the corresponding 
magnitude bin. On the other hand, the STY estimator is ob- 
tained by adopting for ip a parametric function as, for exam- 
ple, a Schechter function (Schechter 1976). For the case un- 
der study we have modified the Schechter function to adapt 
it to the characteristics of our problem: 



<fi{M) oc 10' 



0.4(M* -M){A+1) 



exp(-10 



1.6{M*-M)(A+1) 



) (11) 



M{Mk\e) oc 



Jcp{M\e)dM 



(7) 



where the parameters A and M* are related with the slope 
and with the position of the cut-off of the white dwarf lumi- 
nosity function, respectively. It is worth noticing that this is 
a good characterization of the luminosity function when the 
cut-off is sharp. The real white dwarf luminosity function 
does not show such a sharp cut-off but, instead, a tail ex- 
tending to fainter magnitudes is observationally found (Os- 
walt et al. 1996). This is the reason why this method can- 
not be applied to real samples but to simplified synthetic 
samples in which this tail is not present (see §3 below). In 
both cases, the likelihood can be maximized using standard 
methods. The main drawback of these methods is that they 
can obtain the shape of the luminosity function but do not 
provide the normalization factor. 



2.3 McLximum likelihood estimators based on 
poissonian statistics 

As an alternative to the estimators shown in the previous 

sections, there exist other maximum likelihood estimators 
that build the likelihood using a different approach. Taking 
as a premise that local distribution of objects in some pair 
of variables of the sample has a poissonian distribution, it 
is possible to define a likelihood as a function of the param- 
eter space. The first of such methods, the C~ method, was 
proposed by Lynden-Bell (1971) and was later improved by 
Choloniewski (1986). The Choloniewski method uses simple 
data to define a probabilistic model and then a new likeli- 
hood, by dividing the parameter space (magnitude and par- 
allax in our case) in cells and assuming poissonian statistics 
for each cell. This method allows to estimate both the shape 
of the luminosity function and the total density of objects 
simultaneously. 
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Figure 1. Distance versus absolute magnitude for a sample pop- 
ulation of white dwarfs. Also shown are the limiting magnitude 
and the parameter space, S, adopted for the calculation of the 
luminosity function using the Choloniewski method. See text for 
details. 



We consider a sample with a total number of A'obj ob- 
jects having absolute magnitudes Mi and parallaxes ni, with 
i — 1 , . . . , A^obj . The method is restricted to a solid angle 
f2 = 47r/3. Moreover, the absolute magnitude is assumed 
to fall within M € [Mo, Ma] and the parallax also obeys 
TT € [TrojTTs]. This defines a volume: 

Vt = ^{rl-r'o) (12) 

where ro = I/ttq and vb = l/'Kb- Consequently, the number 
density can be written as 

" = T^ 

and the luminosity function and spatial density of the sam- 
ple are, respectively: 

/ ^{M)dM = n (14) 

J Mo 



(13) 



VtQ{x,y, z) dxdydz = A^obj 



(15) 



where g is the spatial density of the sample of objects. The 
key point of the method is to assume that the number of 
objects in every interval dM dxdydz is given by a poissonian 
random process and that, consequently, the probability of 
finding k objects in each box can be modelled using the 
following expression: 



Pk = exp(-A) 



fc! 



(16) 



If, furthermore, the distribution in absolute magnitude and 
the spatial density distribution are assumed to be indepen- 
dent — which may not be the case if different kinematic 
populations are present in the sample — we have: 

A = -ip(M)q(x, y, z) dM dxdydz (17) 
n 

In order to build the likelihood, the previous expressions 
are discretized by considering the distribution of objects in 
the (M, tt) plane, and binning it into square boxes (see figure 
1). We set AM = Att = A and Mi = Mo -\- iA and nj = 



no + J A with j = 0, . . . , A, and i = 0, . . . , B. We denote the 
number of objects that can be found in the box as Nij 
and so the binned probability can be written as: 

(18) 



PNii = exp(-Aij) 



A 



with: 

Aij = :^ ft A ^ (jj -rj_i) (19) 

being (pi closely related to the luminosity function in the 
given magnitude interval through the expression: 

-1 

dM , (20) 

'Mi J 



Mi -I 



ifii = I <p{M)dM X 
whereas the spatial density is given by 

Qj = / Q{x,y, z)r^ dr dQx 



r dr dil 



(21) 



Finally, we mention that the Choloniewski likelihood 
must be computed considering the selection effects. As the 
sample only provides information for apparent magnitudes 
up to a limiting magnitude m < mum, the value of the num- 
ber of objects in each box is only known for the region below 
the selection line. Thus, we have 



= n n exp(-Ai,) 



Nii 



(22) 



i,jes 



where S stands for grid in the parameter space (see figure 
1). This likelihood can, again, be maximized using standard 
methods. 



3 THE MONTE CARLO SIMULATIONS 

Our Monte Carlo simulator has been thouroughly described 
in previous papers (Garci'a-Berro et al. 1997; Garcia-Berro 
et al. 2004) so here we will only summarize the most im- 
portant inputs. We have used a pseudo-random number 
generator algorithm (James 1990) which provides a uniform 
probability density within the interval (0, 1) and ensures a 
repetition period of >10^*, which is virtually infinite for 
practical simulations. When gaussian probability functions 
are needed we have used the Box-Muller algorithm (Press et 
al. 1986). 

Since we want to test the behaviour of the proposed es- 
timators previously discussed in §2 under different assump- 
tions for the underlying white dwarf population we have 
analyzed a series of different scenarios with controlled stel- 
lar parameters. In a first set of simulations we have adopted 
the most simple prescriptions for the stellar evolutionary in- 
puts. Specifically, we have adopted the most simple cooling 
law (Mostol 1952). Consequently, emission of neutrinos was 
not considered. Crystallization and phase separation were 
also disregarded. Additionally, for all white dwarfs we adopt 
the same cooling sequence, namely that of a typical 0.6 Mq 
white dwarf, independently of its respective mass. Thus, the 
effects of the mass spectrum of white dwarfs arc also com- 
pletely disregarded. The initial-to-final mass relationship 
for white dwarfs and the main sequence lifetime of their 
progenitors adopted here were the analytical expressions of 
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Table 1. Summary of models 



Model 


Cooling 
sequences 


''max 


^-distribution 


1 


Mcstcl (1952) 


250 pc 


Uniform 


2 


Mcstcl (1952) 


1800 pc 


Uniform 


3 


Mestcl (1952) 


1800 pc 


h = 300 pc 


4 


Salaris et al. (2000) 


1800 pc 


Time-dependent h 



Iben & Laughlin (1989). Finally, no bolometric corrections 
were used. Also a very simple galactic model has been used. 
In particular, a standard initial mass function (Scalo 1998) 
and a constant volumetric star formation rate were adopted. 
The velocities have been drawn from normal laws taking into 
account the differential rotation of the disk and the peculiar 
velocity of the Sun with respect to the local standard of rest. 
Since the effect of the spatial distribution of the white dwarf 
population can be important, we have performed two kinds 
of simulations: in the first one a uniform distribution was 
used, whereas for the second one a constant scale height of 
300 pc was assumed. 

In a second stage a set of more realistic model simula- 
tions has also been performed. For this second set of simu- 
lations the cooling sequences of Salaris et al. (2000) which 
incorporate the most accurate physical inputs for the stel- 
lar interior (including neutrinos, crystallization and phase 
separation) and reproduce the blue turn at low luminosities 
(Hansen 1998) have been used. Also, these cooling sequences 
encompass the full range of interest of white dwarf masses, 
so a complete coverage of the effects of the mass spectrum 
of the white dwarf population was taken into account. Be- 
sides, the spatial density distribution is obtained from a scale 
height law (Isern et al. 1995) which varies with time and is 
related to the velocity distributions. 

All the simulations presented here are the average of 
an ensemble of 400 independent realizations. In all the cases 
but in the first one the white dwarf population was modelled 
up to distances of rmax = 1800 pc, in order to avoid the 
effects of the border and were normalized to the local space 
density of white dwarfs within 250 pc, n = 0.5 x lO^** pc^"* 
for Mv < 12.75™'"5 (Liebert, Bergeron & Holberg 2005). 
For model 1 we adopted rmax = 250 pc in order to test the 
effects of a distance-limited sample. Table 1 summarizes the 
main characteristics of each simulation. 



4 RESULTS 
4.1 The 1/Vm, 



method: slope, cut— off and binning 



We start discussing model 1. Since the final goal is to com- 
pute the white dwarf luminosity function a set of restrictions 
is needed for selecting a subset of white dwarfs which, in 
principle, should be representative of the whole white dwarf 
population. We have chosen the following criteria for select- 
ing the final sample: my < 18.5™''^ and n > 0.16" yr~^ as it 
was done in Oswalt et al. (1996). We do not consider white 
dwarfs with very small parallaxes (tt < 0.005"), since these 
are unlikely to belong to a realistic observational sample. 
Additionally, all white dwarfs brighter than My < 13™''^ 




Figure 2. Distance modulus versus absolute magnitude for a 
whole population of white dwarfs within a 1800 pc in a uniform 
distribution. Also shown are the selection criteria in proper mo- 
tion and apparent magnitude. See text for details. 



are included in the sample, regardless of their proper mo- 
tions, since tlie luminosity function of hot white dwarfs has 
been obtained from a catalog of spectroscopically identified 
white dwarfs (Green 1980; Fleming et al. 1986) which is as- 
sumed to be complete (Liebert, Bergeron & Holberg 2005). 
Moreover, all white dwarfs with tangential velocities larger 
than 250 km s~^ were discarded (Liebert et al. 1989) since 
these would be probably classified as halo members. With 
all these inputs the white dwarf luminosity function should 
have a constant slope and a very sharp cut-off, which de- 
pends on the adopted age of the disk. Given the cooling 
law adopted here the slope turns out to be 5/7 and, more- 
over, considering that the cut-off of the observational white 
dwarf luminosity function is located at log (L/Lq) ~ —4.6 
(Liebert et al. 1988) the adopted age of the disk turns out 
to be 13 Gyr. 

In order to illustrate the effects of the previous restric- 
tions in the final sample, in figure 2 we show the distance 
modulus as a function of the absolute V magnitude for the 
whole white dwarf population of a given realization of our 
Monte Carlo simulations. The upper horizontal line corre- 
sponds to the maximum distance to which our synthetic 
population extends (1800 pc). The vertical line corresponds 
to the limiting absolute magnitude, Mum, which is directly 
obtained from the adopted age of the disk and the Mestel 
cooling law used in this set of simulations. The diagonal line 
corresponds to the limiting magnitude imposed by the selec- 
tion criteria previously mentioned. On the other hand, there 
is as well a maximum distance for which a white dwarf could 
be found within the proper motion limit. In particular, an 
object should have a tangential velocity smaller than 250 
km s~^ to be considered as a disk white dwarf, otherwise 
it would be considered a halo member. Note, however, that 
most of the synthetic white dwarfs above this line have tan- 
gential velocities considerably smaller than 250 km s~^ and, 
consequently have proper motions smaller than the proper 
motion cut. We have drawn an horizontal line for this maxi- 
mum distance (r = 330 pc) for which a white dwarf could be 
considered as a disk member. The upper dotted line in this 
diagram corresponds to the distance for which an otherwise 
typical white dwarf with vtan = 30 km s~^ would be included 
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in the final sample of white dwarfs (r = 65 pc). Finally, it is 
worth noticing that the currently available proper motions 
surveys are not sensitive to large proper motions. A repre- 
sentative upper cut in proper motion could be /iu = 2" yr~^. 
The solid bottom horizontal line in Fig. 2 represents the 
corresponding distance for a high velocity white dwarf with 
Vtan = 250 km s~^, representative of the halo population. As 
can be seen in Fig. 2 the effects of the selection criteria are 
dramatical and only an extremely small percentage of the 
whole white dwarf population meets the selection criteria 
and, consequently, are culled for building the white dwarf 
luminosity function. 

Fig. 3a shows 20 independent realizations of the white 
dwarf luminosity function, computed with the 1/Vmax 
method, binned in 2 bins (top panel), 4 bins (middle panel) 
and 8 bins (bottom panel) per decade. The solid lines corre- 
spond to the theoretical expectations previously described. 
That is, a straight line with constant slope and a very sharp 
cut-off at the observed position. Each sample typically con- 
tains about 300 white dwarfs, the size of the observational 
sample. We recall that, by construction, our samples are 
complete. As can be seen, there is a considerable spread 
about the theoretical expectations and, moreover, the white 
dwarf luminosity function is underestimated at moderately 
high luminosities. On the other hand. Fig. 3 b shows the av- 
erage and standard deviation of the 400 realizations of the 
white dwarf luminosity function. Again it is clearly visible 
that the 1/Vmax method considerably underestimates the 
density of white dwarfs with moderately high luminosities. 
In fact, the 1/Vmax method only recovers the right slope 
for luminosities smaller than log(L/LQ) ~ —2.2. Moreover, 
the position of the cut-off also depends on how the data 
is binned, its position being more accurate for finer bin- 



ning, and it is always located at larger luminosities, a direct 
consequence of the binning procedure, as already found by 
Garci'a-Berro et al (1999). 

In order to quantify the previous statements Fig. 4a 
shows the frequency distribution of slopes for the 400 inde- 
pendent realizations of the white dwarf luminosity function. 
The vertical solid line corresponds to the theoretical value of 
the slope of the white dwarf luminosity function (5/7). Ob- 
viously all the realizations severely overestimate the slope 
and, moreover, there is a considerable dispersion. On the 
other hand, in Fig. Qf) the distribution of simulated cut-offs 
is shown. Clearly, the finer the binning the more accurate is 
the determination of the cut-off. However, the dispersion is 
relatively small. That is, this is a systematic bias, which can 
be accounted for and ultimately corrected. 

Another important information that can be readily ob- 
tained from our simulations is how to compare the observa- 
tional procedure to assign error bars — basically assuming 
poissonian statistics for each bin (Liebert et al. 1988) — 
with the computed standard deviations, Alogi^std- To be 
more specific the contribution of each star to the total er- 
ror budget in its luminosity bin, A log i/^obs , is conservatively 
estimated to be the same amount that contributes to the re- 
sulting density; the partial contributions of each star in the 
bin are squared and then added, the final error is the square 
root of this value. Table 2 shows the result of such a com- 
parison. As can be seen, we have found that, in general, the 
standard procedure to assign error bars severely underesti- 
mates the observational error bars, especially for the brigth- 
est luminosity bins of the luminosity function. Particularly, 
for these luminosity bins the error bars are underestimated 
by a factor of roughly « 10, whereas for the two (more 
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Figure 4. a) Distribution of slopes for the white dwarf luminosity functions of Fig 3. b) Distribution of cut-offs for the white dwarf 
luminosity functions of Fig. 3. 



Table 2. A comparison of the error bars of the white dwarf lumi- 
nosity function computed using the 1/Vmax method, Alogi^obsi 
with the standard observational procedure to assign error bars 
with the intrinsic statistical deviation of the 400 realizations for 
the white dwarf population of model 1, Alogipstd. The white 
dwarf luminosity function has been obtained by binning the syn- 
thetic data in four luminosity bins per decade. 



-log(L/LQ) 


A log Ifobs 


A log (fistd 


1.625 


-7.155 


-5.796 


1.875 


-7.222 


-5.535 


2.125 


-7.097 


-5.207 


2.375 


-5.886 


-5.016 


2.62-5 


-5.824 


-4.8,38 


2.875 


-6.000 


-4.689 


3.125 


-5.699 


-4.405 


3.375 


-5.398 


-4.291 


3.625 


-5.046 


-4.165 


3.875 


-5.155 


-3.925 


4.125 


-4.823 


-4.064 


4.375 


-4.921 


-5.063 



populated) dimmest luminosity bins the error bar and the 
inherent statistical deviation are very similar. 



4.2 Looking for an alternative: other estimators 

As we have shown, the 1/Vmax method does not provide sat- 
isfactory answers with regard to the slope of the white dwarf 
luminosity function and the position of the cut-off for the 
white dwarf populations of model 1. Thus, some alternatives 
must be explored. As previously mentioned, there exist sev- 



eral such alternatives. For the sake of conciseness here we 
will discuss only two of them: the STY method (Sandage et 
al. 1979) and the Choloniewski method (Choloniewski 1986). 
The SWML method (Efstathiou et al. 1988) gives results 
which are very similar to the STY method and, thus, we 
will not describe them in detail for the moment. AH three 
methods are maximum likelihood methods and have been 
consistently used to estimate galaxy luminosity functions. 
And this is perhaps their main drawback since they have 
not been devised to correct for the bias in proper motion. 
This is a characteristic of the current white dwarf samples, 
and the 1/Vmax method does correct it. However, as it will 
be shown below this is not a severe problem, at least for the 
STY method. We must recall that the STY method pro- 
vides the shape of the luminosity function but not its nor- 
mahzation (namely, the true density of objects), whereas 
the Choloniewski method provides both the shape of the 
luminosity function and the density of objects. The selec- 
tion criteria used in this set of simulations arc exactly the 
same used previously for deriving the white dwarf luminosity 
function using the 1/Vmax method. 

Fig. 5a shows a comparison of the results obtained for 
our model 1 using the different methods discussed here. The 
bottom panel is the same already shown in the middle panel 
of Fig. 3 b. The middle panel shows the results obtained us- 
ing the STY method, and as can be seen, the STY method 
recovers very well the correct slope. Finally, the top panel 
shows the results obtained using the Choloniewski method. 
Clearly, this method underestimates the slope at high lu- 
minosities. The reason for this is quite clear: the statistics 
of the brightest luminosity bins are not poissonian. In fact, 
very few white dwarfs populate these bright luminosity bins 
(note that the vertical scales in Fig. 3 and 5 are logarithmic) . 
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Figure 5. a) A compajrison of the white dwarf luminosity functions obtained using alternative methods and the 1/Vmax method, b) 
Distribution of cut-offs for the different realizations of the white dwarf population using different methods. 



Consequently, the results obtained using the Choloniewski 
method for the brightest luminosity bins are not correct. 
However this method turns out to be very useful since we 
obtain the correct density of white dwaxfs. In all three cases 
the error bars are similar. On the other hand, Fig. 5b shows 
the position of the cut -off for all three methods. The STY 
method provides better results than the 1/Vmax method but, 
undoubtely, the Choloniewski method performs the best, al- 
though with a larger variance than the STY method. Finally, 
the statistical error bars for all three methods are rather sim- 
ilar, being somehow smaller for the Choloniewski method. 



4.3 Extending the sample to larger distances 

One possible reason for the systematic bias found when us- 
ing the 1 /Vmax method to obtain the white dwarf luminos- 
ity function could be due to the fact that bright objects can 
be found at distances considerably larger than 250 pc, the 
maximum distance within which we have distributed syn- 
thetic white dwarfs in model 1. For instance, an object with 
\og{L / Lq) = —1.0 will be within our apparent magnitude 
selection criterion oven if it is located at distances as far 
as 1800 pc. For this reason we have applied the different 
luminosity function estimators to a sample with a larger 
maximum distance, in particular to a sample of synthetic 
white dwarfs distributed within a sphere of radius 1800 pc 
(model 2 in Table 1). Since in this sample the effects of the 
scale height should be important we have carried out an 
additional simulation in which the synthetic white dwarfs, 
instead of being distributed according to an uniform density 
law, have been distributed according to an exponentially de- 
crasing density law. This model will be discussed in section 

4.4 below. In both cases the rest of the galactic and stellar 
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Figure 6. Same as Fig. 5a for the white dwarf population of 
model 2. 



evolutionary inputs remain unchanged. Note, however, that 
in the sample of spectroscopically identified white dwarfs of 
McCook & Sion (1999) — the primary observational source 
from which the white dwarf luminosity function is built — 
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Figure 7. Distance modulus versus magnitude for a typical re- 
alization of the white dwarf population of models 1 and 2. The 
diagonal lines represent the adopted bins of the white dwarf lu- 
minosity funcion. The top panel shows the results for a sphere 
of 1800 pc with a uniform density law (model 2), whereas the 
bottom panel shows the results for a sphere of 250 pc, also with 
a uniform density distribution. 



the most distant white dwarf with a reUable parallax deter- 
mination is located at ~ 250 pc. 

In figure 6 the luminosity functions for the white dwarf 
population of model 2 are shown for the three estimators 
previously discussed. The results shown here are the ensem- 
ble average and standard deviation of a set of 20 realiza- 
tions. As can bo soon, now the 1/Vmax estimator correctly 
matches the theoretical expectations for the slope of the 
white dwarf luminosity, although with a considerably large 
variance for the brightest luminosity bins. The reason for 
this behavior will be discussed below, with the help of figure 
7. Also, it is interesting to note that for this set of sim- 
ulations the STY estimator also yields reasonable results, 
whereas the Choloniewski estimator largely underestimates 
the white dwarf density for the luminosity bins near the 
maximum of the white dwarf luminosity function. Finally, 
the three estimators obtain the same cut-offs previously ob- 
tained for model 1, as it should be expected, given that the 
population of faint white dwarfs is drawn from distances 
smaller than 300 pc. 

In Fig. 7 the distribution of the distance modulus as a 
function of the magnitude for the synthetic white dwarf pop- 
ulations of models 1 (bottom panel) and 2 (top panel) are 
shown. We have also represented the diagonal lines corre- 
sponding to the adopted bins of the white dwarf luminosity 
function, from log(L/LQ) = —1.0 (top line) to —5.0 (bot- 
tom line). The horizontal line corresponds to r = 250 pc. 
As can be seen, for model 2 (the sample distributed within 
1800 pc) a sizeable number of intrinsically bright white 
dwarfs (at largo distances) meet the selection criteria and, 
consequently, contribute to the white dwarf luminosity func- 
tion. In the sample obtained from model 1 these intrinsically 
bright white dwarfs are missing, and this is the reason why 
we obtain a biased white dwarf luminosity function. It is nev- 
ertheless worth mentioning three important points. Firstly, 
most of the spectroscopically identified white dwarfs in the 
catalog of McCook & Sion (1999) have distances smaller 
than 250 pc, as can be seen in top panel of figure 8. In fact. 
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Figure 8. Same as figure 7 for model 3 and for the sample of spec- 
troscopically identified white dwarfs of McCook &; Sion (1999). 



only three white dwarfs in this catalog have distances larger 
than 250 pc and, conscqucutly, the slope of the bright branch 
of the white dwarf luminosity function obtained using this 
catalog as the primary source of observational data should 
be, in principle, questioned. Secondly, the bright branch 
of the white dwarf luminosity function depends primarily 
on the relative strengths of neutrino leakage and radiative 
losses. Hence, a robust determination of the slope of the 
bright branch of the white dwarf luminosity function turns 
out to be important for deriving very useful constraints on 
the physics of cooling white dwarfs. And, finally, and most 
importantly, even in the case of a distance limited sample, 
such as that of model 1, the maximum likelihood estimar 
tors and, specifically the STY method and the Choloniewski 
method, detect the deficit of intrinsically bright white dwarfs 
and, moreover, they are able to retrieve the correct slope. 
Thus, these estimators are much more robust and provide 
more reliable white dwarf luminosity functions. 



4.4 The effects of the scale height 

We have also tested which would the depedence of the white 
dwarf luminosity function on an exponentially decreasing 
scale height law. Hence, instead of assuming an uniform 
distribution of white dwarfs within the computational vol- 
ume, as it has been done so far for models 1 and 2, we have 
adopted a constant scale height of 300 pc (model 3 in Table 
1), and we have distributed our synthetic white dwarfs ac- 
cordingly within a volume of radius 1800 pc. The results for a 
typical realization of our Monte Carlo simulations are shown 
in the bottom panel of figure 8, where the distribution of dis- 
tance modulus as a function of the apparent magnitude is 
shown for this synthetic population. As can be seen, for the 
brightest luminosity bin we obtain now very few synthetic 
stars, as it should be expected, given that the number den- 
sity of white dwarfs at large distances is heavily suppressed 
by the adopted exponentially decreasing scale height law. 
In fact we now obtain more or less the same number of 
stars in both the synthetic sample and in the observational 
sample of McCook & Sion (1999). Specifically, for model 3 
we obtain on average three synthetic white dwarfs in the 
brightest luminosity bin, whereas in the observational sam- 
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pie of McCook & Sion (1999) two white dwarfs are found 
in this luminosity bin. These numbers should be compared 
with the corresponding one for the population of synthetic 
white dwarfs obtained for model 2, which, on average, is 13 
synthetic stars, significantly larger than the previous ones. 
However, it is worth noticing as well that the effects of the 
scale height and of the completeness of the sample under 
study — especially at large distances — are difficult to dis- 
entangle, at least for the observational sample of McCook 
& Sion (1999). Clearly, a much better observational catalog, 
complete up to distances considerably larger than the scale 
height of Galactic disk, should be needed to this regard. For 
the moment being, the possibility of obtaining the value of 
the scale height from a sample of tipicaly 300 white dwarfs 
remains very remote. 

4.5 A realistic model 

As a final test we have applied the several estimators dis- 
cussed here to a more realistic model of the white dwarf 
population, which we denote as model 4 (see Table 1). We 
recall that for model 4 we have used the cooling sequences 
of Salaris et al. (2000), which encompass the full range of 
masses of interest, as opposed to what has been done up to 
now, where the cooling rate of any white dwarf, regardless 
of its mass, was obtained from a single cooling sequence of a 
0.6 Mq white dwarf. Moreover, these cooling sequences in- 
corporate the effects of neutrinos, crystallization and phase 
separation. Consequently the slope of the white dwarf lumi- 
nosity function is no longer constant but, instead, reflects 
the relative speed of cooling at a given luminosity. In partic- 
ular, for those luminosities where neutrino cooling is domi- 
nant the cooling rate is larger and, consequently the slope of 
the white dwarf luminosity function turns out to be steeper, 
yielding less white dwarfs for these luminosity bins when 
compared to the fiducial luminosity function used up to 
now. Conversely, for those luminosities where crystallization 
and phase separation are the relevant physical processes, the 
cooling speed is smaller (the release of crystallization latent 
heat and the gravitational energy released by phase separa- 
tion must be radiated away) and, consequently, the slope of 
the luminosity function is also steeper, yielding in this case 
more white dwarfs for these luminosity bins than the fidu- 
cial luminosity function obtained from Mestel's law (since 
they pile-up at these luminosities due to a reduced cooling 
rate). For this reason, the expression of Eq. (11) for the 
STY estimator is no longer valid and, consequently, instead 
of using the STY estimator for this set of simulations we 
adopt the SWML method, which provides a more appro- 
priate computational approach. We also note that in this 
case we have adopted our full model of Galactic evolution, 
as described in detail in Garci'a-Berro et al. (1999). Within 
this model the adopted scale-height depends on time — be- 
ing larger for past epochs — and, consequently, since the 
adopted star formation rate in the local column has been 
adopted to be constant the volumetric star formation rate is 
no longer constant. Moreover, the velocity dispersions also 
depend on time and, thus, the distributions of velocities are 
not perfectly gaussian as it was the case for models 1 to 3. 
As a matter of fact our Galactic evolutionary model natu- 
rally incorporates the thin and the thick disk populations 
— see Torres et al. (2002). However, the faint end of the 



disk white dwarf luminosity function is generally assumed 
to be contaminated by a yet not well known fraction of halo 
white dwarfs (Reid 2005). Indeed, although the peak of the 
halo white dwarf luminosity function is located at a lumi- 
nosity considerably fainter than that of the cut-off of the 
disk white dwarf luminosity function (Isern et al. 1998) some 
halo white dwarfs may be present in faintest luminosity bins. 
This is the reason why we apply a very strict velocity cut of 
250 km s~^. While it is true that this simple procedure does 
not completely remove high velocity populations it is also 
true that the results obtained with model 4 represent a step 
in the right direction. Finally we point out that in order to 
keep consistency with the simulations previously described 
we have adopted the same age of the disk. Since the cooling 
sequences of model 4 incorporate the effects of crystalliza- 
tion and phase separation, which introduce a sizeable delay 
in the cooling times, the cut-off in the white dwarf luminos- 
ity functions moves to fainter luminosities accordingly. 

At this point of the discussion of our results it is im- 
portant to realize that up to now we have always had a 
"template" white dwaxf luminosity function to which we 
could compare. This template was the very simple lumi- 
nosity function already shown in Figs. 3, 5 and 6. Given 
the stellar and galactic inputs adopted for model 4, a white 
dwaxf luminosity function with a perfectly constant slope 
and a sharp cut-off is a poor characterization of the theo- 
retical expectations. However, we can easily obtain a tem- 
plate white dwarf luminosity function in the following way. 
We recall that, by construction, our samples are complete, 
although we only select about 300 white dwarfs using the se- 
lection criteria already discussed before. However, our sim- 
ulations do provide the whole population of white dwarfs. 
Hence, we can obtain the reaJ luminosity function by simply 
counting white dwarfs in the computational volume. This 
is done for all realizations and then we obtain the average. 
The result is depicted as a solid line in Fig. 9, where we also 
show the results obtained using the Choloniewski method 
(upper panel), the SWML method (middle panel) and the 
1/Vmax method (bottom panel), with their computed error 
bars. As can be seen the cut-off of the white dwarf lumi- 
nosity function has moved to fainter luminosities, its precise 
location being now log(Z//Z/Q) ~ —4.8, a direct consequence 
of crystallization and phase separation. 

Fig. 9 clearly shows that the performance of the 1/Vmax 
method is superb, since this method nicely flts both the 
shape of the white dwarf luminosity function and the posi- 
tion of cut-off. The SWML method (middle panel of Fig. 
9) also flts pretty well the shape of the white dwarf lumi- 
nosity function, but the last two (faint) luminosity bins are 
poorly determined. Consequently, the determination of the 
cut-off of the white dwarf luminosity function is subject 
to a large variance, and individual simulations can yield 
very different results for the age of the disk. Finally, the 
Choloniewsky method (top panel of Fig. 9) clearly underes- 
timates the number of faint white dwarfs (the peak in the 
white dwarf luminosity function) and does not reproduce 
the real cut-off. AH in all, for a sample of about 300 white 
dwarfs, and when all the observational biases are correctly 
taken into account the 1/Vmax performs best. 

Also, some computational details are worth mention- 
ing. The first one is that the computational load of the two 
maximum likelihood methods is much larger than that of 
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Figure 9. White dwarf luminosity function for model 4 using the 

different estimators under study (dots), our template luminosity 
function is shown as a solid line. 
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Figure 10. Same as Fig. 9 for a sample of 2000 white dwarfs, the 
sample size of the SDSS. 



the 1/Vmax method. This does not pose a severe problem 
when samples with a small number of objects axe analyzed 
but it is a point to be considered when samples containing 
a large number of white dwarfs, like that of the SDSS which 
will be the object of §4.6, are studied. The second impor- 
tant remark is that for a sample size of 300 white dwarfs 
the convergence of the two maximum likelihood methods is 
slow, a consequence of the minimum being very shallow. 



4.6 The future: the SDSS 

Very recently, a sample of white dwarfs selected from the 
Sloan Digital Sky Survey Data Release 3 (SDSS DR3) com- 
bined with improved proper motions from the USNO-B has 
derived a preliminary (although very much improved) white 
dwarf luminosity function based on 6000 stars (Harris et al. 
2006). We emphasize at this point that we do not aim to per- 
form a full analysis of the sample of Harris et al. (2006), but 
a preliminary assessment of it. A detailed analysis of this 
sample is out of the scope of this paper and we postpone 
it for a forthcoming publication. The white dwarf luminos- 
ity function of Harris et al. (2006) has been built using the 
following selection criteria. The survey area of the SDSS is 
mostly centered around the North Galactic Cap and cov- 

2 

ers an area of 5282 . For our Monte Carlo simulations we 
have adopted the precise geometry of the SDSS, an ellipti- 
cal region centered at a = 12'' 20™ , 6 = +32.8°, whose 
minor axis is the meridian at that right ascension, with ex- 
tent ±55° in declination. The major axis is the great circle 
perpendicular to that, and the extent is ±65°; it extends 
from about 7^ e"^'"" to about 17''34""'°. Prom the original 



sample of 6000 stars Harris et al. (2006) have only selected 
stars with /i > 20 mas yr^^ and, thus, we disregard all white 
dwarfs with proper motions smaller than this value. Addi- 
tionally, Harris et al. (2006) use the reduced proper motion 
Hg = g + 51ogAt + 5 = Mg + 5 log Kan - 3.379, where g 
is the SDSS magnitude, to discriminate between main se- 
quence stars and white dwarfs, since the latter are typically 
5-7 magnitudes less luminous than subdwarfs of the same 
color. Moreover, they require that all white dwarfs must have 
Vtan > 30 km/s to enter in the final sample, and this is what 
we adopt. An additional criterion is that all white dwarfs 
should have 15.0 < g < 19.5. We have selected only white 
dwarfs with 15.0 < my < 19.5. The final size of the sample 
used to built the white dwarf luminosity function is of about 
- 2000 stars. 

With all these restrictions we have computed the white 
dwarf luminosity function with the inputs of model 4. The 
results are shown in Fig. 10. As it has been done so far, we 
show the white dwarf luminosity function computed with the 
Choloniewski method (top panel) , the SWML method (mid- 
dle panel) and the 1 /Vmax method (bottom panel) . Clearly, 
both the Choloniewsky method and the 1/Vmax method per- 
form very well, whereas the SWML method misses the max- 
imum and the cut-off of the white dwarf luminosity function 
and, moreover, the variance for the brightest luminosity bins 
is much larger than those of the other two methods. For the 
Choloniewski method the last luminosity bin does not show 
up, but it should be taken into account that the the vari- 
ance of the last bin of the 1/Vmax method is very large. 
One comment is in order regarding this last method. As can 
be seen in Fig. 10, the 1/Vmax method underestimates the 
white dwarf luminosity function for the brightest luminos- 
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Figure 11. Same as Fig. 10 but assuming now an incomplete- 
ness of the input catalog of 20% (filled symbols) and 40% (open 
symbols). The open symbols have been moved by Alog(L/L0) = 
—0.08 for the sake of clarity. 

ity bins. This is a consequence of the adopted galactic in- 
puts for our white dwarf population and, more specifically, 
of the adopted scale height. Since we are using the original 
1/Vmax method, whithout correcting for the scale height, 
the number of white dwarfs per unit volume and magnitude 
interval is underestimated for the brightest luminosity bins, 
where the survey extends to relatively largo distances. On 
the other hand, the Choloniewski method overestimates the 
white dwarf density for these luminosity bins. Note however, 
that for the intermediate luminosity bins the 1 /Vmax method 
matches very well the shape of the white dwarf luminosity 
function. All in all, except for the brightest luminosity bins, 
the 1/Vmax method provides a very good characterization of 
the white dwarf luminosity function. Finally, and contrary 
to what was found in §4.1 for a sample of 300 white dwarfs, 
the observational procedure for assigning the error bars to 
the white dwarf luminosity function is fair for a sample of 
2000 white dwarfs. 



4.7 The incompleteness of the sample 

Another important concern is how the incompleteness of the 
sample affects the shape and the location of the cut-off of 
the retrieved white dwarf luminosity function, and how ro- 
bust are the different methods when a sizeable fraction of 
the input sample is discarded. This is precisely the goal of 
this section. In order to assess these effects wo have first ran- 
domly eliminated from the final input sample discussed in 
§4.6 (that is the sample simulating the white dwarf luminos- 
ity function obtained from the SDSS DR3) 20% and 40% of 
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Figure 12. Differences of the resulting white dwarf luminosity 
function, Alog</), when incompletenesses of 20% (solid squares), 
40% (open squares) and a linealy decreasing completeness (solid 
triangles) are assumed, with respect to the white dwarf luminosity 
function obtained when the full input sample is used. See text for 
details. 

the white dwarfs which pass all the selection criteria, inde- 
pendently of their magnitude, proper motion, or tangential 
velocity. The results arc shown in Fig. 11. 

As can be seen in the top panel of Fig. 11, the white 
dwarf luminosity functions obtained using the Choloniewski 
method do not differ considerably from those previously 
studied in §4.6 and, consequently, this method is extremely 
robust against possible incompletenesses of the input sam- 
ple, even under the radical assumption that about 40% of 
the white dwarfs in the input sample are discarded in the se- 
lection process or, simply, missed whatever the cause could 
be. For the case in which the SWML method is used we 
stress that this method has the shortcomings already com- 
mented before: firstly, it only recovers the shape of the lu- 
minosity function but not the total density of white dwarfs, 
and, secondly, it misses the faint end of the white dwarf lu- 
minosity function. However, there are not big differences in 
the recovered shape of the white dwarf luminosity function, 
even for incompletenesses of the order of 40%. This is not 
the case of the 1 /Vmax method which, as can be seen in the 
bottom panel of Fig. 11, largely underestimates the result- 
ing white dwarf density for almost all the luminosity bins. 
Note, however, that in this case the luminosity of the cut-off 
is correctly retrieved, independently of the adopted incom- 
pleteness. Hence, our results show that the Choloniewski 
method is much more stable than the 1/Vmax method, oven 
under extreme assumptions about the completeness of the 
input sample used to build the white dwarf luminosity func- 
tion. On the other hand, for the case of the 1/Vmax method 
the size of the error bars is more or less the same than in 
the case in which the sample was complete. This is not the 
case for the brightest luminosity bins when the Choloniewski 
method is used. 

In a second step we have adopted a different strategy. 



The white dwarf luminosity function. I. Statistical errors and alternatives 13 



Instead of discarding a given percentage of white dwarfs re- 
gardless of their properties we have assumed than the in- 
put sample is complete for apparent magnitudes mv = 15.0 
and that the completeness decreases lineraly to 60% for 
mv = 19.5. Fig. 12 shows the results of applying this proce- 
dure to the input sample. For this figure we have preferred 
to show the differences A log ip = log ip' — log ip of the re- 
sulting luminosity function, log ip' , with respect to the white 
dwarf luminosity function obtained using the full input sam- 
ple, log p, in order to better visualize the results. The solid 
squares are the differences when a completeness of 77 = 80% 
is assumed, the open squares are the data for a completeness 
of only 60% and the triangles represent the results obtained 
when a linarly decreasing completeness is adopted. Fig. 12 
shows that the 1/Vmax method underestimates the white 
dwarf luminosity function for the vast majority of the lu- 
minosity bins for all three cases, whereas the Choloniewski 
method is quite robust and, except for the brightest luminos- 
ity bins, is rather insensitive to the completeness of the input 
sample. Hence, and from this point of view the Choloniewski 
method is clearly superior to the 1/Vmax method. 



5 CONCLUSIONS 

We have performed a study of the statistical reliability of 
the white dwarf luminosity function using different esti- 
mators. These include the classical 1/Vmax method, and 
two parametric maximum-likelihood estimators, namely the 
Choloniewski method and the SWML or the STY method, 
depending on the adopted cooling sequences. In a first stage, 
for all three estimators the input sample was drawn from a 
controlled sample for which we adopted the most simple 
cooling law (Mestel 1952) and very schematic galactic in- 
puts. This was done in order to study the real behavior of 
the estimators and to isolate their respectives advantages 
and drawbacks. Nevertheless, for these numerical experi- 
ments the observational selection criteria were fully taken 
into account. We have found that for a small sample size 
the 1 /Vmax method provides a poor characterization of the 
bright end of the white dwarf luminosity function if the sam- 
ple selection procedure is not done carefully. Specifically, this 
method produces an artificial deficit of white dwarfs at mod- 
erately high luminosities when the sample does not contain 
intrinsically bright white dwarfs located at relatively large 
distances. This is a direct consequence of the scarcity of 
intrinsically bright white dwarfs which, in turn, is a conse- 
quence of the very short evolutionary time-scales of these 
white dwarfs. We have, furthermore, shown that this is pos- 
sibly the case of the catalog of spectroscopically identified 
white dwarfs of McCook & Sion (1999), for which very few 
intrinsically bright white dwarfs are present. Moreover, we 
have also demonstrated that for a sample size of 300 stars, 
the 1/Vmax method overestimates the position of the drop- 
off of the white dwarf luminosity function. This is a conse- 
quence of the small number of objects in the input sample 
which, in turn, forces a coarse binning. We have further dis- 
cussed the effect of the adopted scale height law and we 
have found that for a sample size of 300 stars its effect can- 
not be disentangled from the effects of the sample selection 
procedure. Additionally, we have also shown that the obser- 
vational procedure to assign error bars is too optimistic for 



small sample sizes, with realistic error bars being typically 
10 times larger for a typical sample size of 300 objects. 

We have explored two alternatives, the STY method 
and the Choloniewski method. Both methods have been 
widely used to build galaxy luminosity functions with sat- 
isfactory results, and we have found that for the case of 
small sample sizes they perform considerably better than 
the 1/Vmax method, even if none of the two methods takes 
into account the bias of proper motion selected samples. In 
particular, the STY method performs best at recovering the 
slope of the luminosity function whereas the Choloniewski 
method recovers best the position of the cut-off. However, 
the STY method does not provide the true density of white 
dwarfs, whereas the Choloniewski method does. 

We have also applied the two maximum likelihodd 
methods and the 1 /Vmax method to a sample of 300 white 
dwarfs obtained using realistic stellar and galactic inputs. 
In this case, instead of using the STY method the SWML 
method was used, given that the slope of the increasing 
branch of the white dwaxf luminosity function is no longer 
constant. We have found that all three methods present large 
variances for the brightest luminosity bins, that the SWML 
method and the 1 /Vmax method retrieve the correct location 
of the cut-off and that the Choloniewski method underesti- 
mates the number of faint white dwarfs, resulting in a bad 
characterization of the maximum and of the cut-off of the 
white dwarf luminosity function. 

Finally, we have also applied these three methods to a 
sample of 2000 white dwarfs, which is representative of the 
sample used to build the white dwarf luminosity function 
from the SDSS DR3 (Harris et al. 2006). This input sam- 
ple was obtained using up-to-date cooling sequences, real- 
istic galactic inputs and an accurate sample selection pro- 
cedure, following very precisely the prescriptions used for 
drawing the final sample of white dwarfs of the SDSS DR3. 
We have found that the performances of the Choloniewski 
method and of the 1/Vmax method are very similar, pro- 
viding with reasonable accuracy both the detailed shape of 
the white dwarf luminosity function and the location of the 
cut-off. Consequently, in principle both methods could be 
used in a real case, yielding similar results. On the other 
hand, the SWML method does not recover neither the cor- 
rect shape of the luminosity function nor the position of the 
cut-off and, consequently, should not be used for a real sam- 
ple. We have also demonstrated that the effects of the scale 
height law are non-negligible for both the Choloniewski and 
the 1/Vmax method. Particularly, this last method under- 
stimates the white dwarf density for the brightest luminos- 
ity bins, whereas the Choloniewski method overestimates it. 
For this last input sample we have also analyzed the effects 
of the incompleteness, finding that only the Choloniewski 
method is robust when the possible incompleteness of the 
sample is taken into account, retrieving the correct total 
density of white dwarfs even for severe incompletenesses of 
the input sample. In particular, the 1/Vmax method severely 
underestimates the total number density of white dwarfs for 
sample sizes of the order of 2000 stars when an incomplete- 
ness of 20% is adopted, whilst the Choloniewski method 
does not, being thus much more robust than the classical 
1/Vmax method. However, the 1/Vmax method nicely recov- 
ers the position of the cut-off of the white dwaxf luminos- 
ity function, whereas the Choloniewski method does not. 
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In summary, when the input sample has a sizeable number 
of objects a combination of both the Choloniewski and the 
1/Vmax method provides reliable determinations of the white 
dwarf luminosity function. Other estimators, like the SWML 
method, axe not recommended whatsoever given that, firstly, 
they do not provide the true density of white dwarfs but only 
the shape of the luminosity function and, secondly, they do 
not have a performance better than the other two methods. 
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